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We study the phase diagram of the half-filled one-dimensional extended Hubbard model at weak 
coupling using a novel functional renormalization group (FRG) approach. The FRG method includes 
in a systematic manner the effects of the scattering processes involving electrons away from the Fermi 
points. Our results confirm the existence of a finite region of bond charge density wave (BCDW), 
also known as a "bond order wave" (BOW), near U = 2V and clarify why earlier g-ology calculations 
have not found this phase. We argue that this is an example in which formally irrelevant corrections 
change the topology of the phase diagram. Whenever marginal terms lead to an accidental symmetry, 
this generalized FRG method may be crucial to characterize the phase diagram accurately. 
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The one-dimensional extended Hubbard model (EHM) 
has been studied extensively for many years, both be- 
cause of its rich phase diagram ,l] and because of its 
possible applications to quasi-lD organic crystals |2| and 
conducting polymers Q . Despite this long history, a con- 
troversy has recently arisen concerning the possible exis- 
tence of a bond order charge density wave (BCDW, also 
called a "bond order wave" (BOW)) phase separating the 
well-known spin density wave (SDW) and charge den- 
sity wave (CDW) phases of the EHM at half-filling. This 
phase has been suggested by Nakamura and supported 
by quantum Monte Carlo (QMC) and more recently 
density matrix renormalization group (DMRG) 0] calcu- 
lations. However, this phase was not obtained in earlier 
numerical and analytical work 0, H, El EH EH EH E3 • 
In particular, this phase is absent in standard one-loop 
g-ology 7, 8] and bosonization calculations. This dis- 
agreement poses a serious question, which we elucidate 
here. We reconcile the recent numerical results 0, 
with g-ology by introducing a functional generalization 
of the standard g-ology formalism. Our functional renor- 
malization group (FRG) method offers a consistent and 
well-controlled approximation that predicts a finite re- 
gion in parameter space in which the BCDW phase is 
spontaneously formed for the half-filled EHM at weak 
coupling. Our results thus go beyond the important ear- 
lier study of Tsuchiizu and Furusaki (14] , who were able 
to obtain the BCDW phase with a standard RG using ad 
hoc approximations. 

The Hamiltonian of the EHM is given by 

H = -t^(cl +la Ci i<T + H.c,) + U'J2ni i ' t n itl 

+V^nin i+1 - ^L^rii, (1) 

i i 
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where t, U, and V are the nearest-neighbor hopping, 
on-site interaction, and the nearest-neighbor interaction, 
respectively. Here we study the EHM at half-filling 
(/j = 0). It is well established that for repulsive in- 
teractions (U, V > 0), the system is in a CDW phase 
for large values of V/U and in a SDW phase for small 
V/U. Weak-coupling RG studies 0, ||| find the bound- 
ary between these two phases to be at U = 2V. Early 
strong-coupling numerical studies and higher-order 
perturbation theory [lfj have found the phase boundary 
to be slightly shifted away from the U = 2V line, with 
a larger SDW phase. Stochastic series expansion QMC 
studies || found that the BCDW phase exists in a finite 
region around the line U = 2V and that it ceases to ex- 
ist when the interaction exceeds a critical value. There 
are disagreements between the published DMRG results. 
An earlier result showed that the BCDW exists only 
precisely at the CDW/SDW phase boundary at interme- 
diate couplings. A more recent DMRG calculation lfj 
obtained the same phase diagram as the QMC study |f| . 

The standard RG-" g-ology" - has proven to be a pow- 
erful method for studying low-energy properties of inter- 
acting one-dimensional systems at weak-coupling. By in- 
tegrating out high-energy modes, one obtains flow equa- 
tions for the marginal couplings such as the two body in- 
teraction vertices. These interaction vertices are in prin- 
ciple functions of three momenta, which can take any 
value within the Brillouin zone and correspond to the 
momenta of the incoming and outgoing electrons. The 
fourth momentum is determined by momentum conser- 
vation. In standard g-ology, the interaction processes are 
classified according to the branch label (right or left) of 
the electrons involved. All further dependence on the 
momenta, i. e., the dependence on the magnitude of 
the momenta, is neglected, since only the dependence on 
the direction of the momenta is marginal. The radial 
dependence is irrelevant according to scaling and power- 
counting arguments |l8j . It is important to notice that 
irrelevant operators renormalizc to zero as the RG pro- 



2 



ceeds but may not be small in the beginning of the flow. 
This is the key issue here, and we will return to it later. 

In g-ology, interaction processes are classified into 
backward scattering (<?i), forward scattering involving 
electrons from two branches (g 2 ) and from the same 
branch (34), and Umklapp process (53)- As the scat- 
tering between electrons with the same spin can be ob- 
tained from scattering between electrons with different 
spins we shall ignore all spin indices, leaving it un- 
derstood that all processes are between electrons with 
different spins. For the EHM, the bare values of the cou- 
plings are g\ = g 3 = U — 2V and g 2 — g± = U + 2V. 

Exactly at U — 2V, both gi and g 3 are equal to 
zero, and they remain zero under the RG flow, result- 
ing in a massless theory for both the spin and charge 
sectors. This is the underlying reason that conventional 
weak-coupling calculations (both g-ology and bosoniza- 
tion) find a direct transition between the SDW and the 
CDW phase exactly at U = 2V, where both gaps van- 
ish simultaneously. An important insight was provided 
by Nakamura £| (and further explored by Tsuchiizu and 
Furusaki 0]), wn0 observed that there is no symmetry 
principle that enforces g\ and 173 to vanish simultaneously 
and that higher-order corrections may lift this degeneracy 
and thereby change the topology of the phase diagram 
151. They then adopted an idea from Penc and Mila 
H| and applied the following two-step procedure: for the 
high-energy part of the band (A > A cuto ff), second-order 
perturbation theory is performed to find corrections for 
the couplings g; these values are then used as the initial 
conditions for the RG procedure, which is performed for 
the low energy part of the band (A < A cuto ff). This 
is sufficient to generate a finite region of BCDW phase. 
Clearly, this procedure is ad hoc and relies on an arbi- 
trary choice for A cuto ff (which in is chosen to be 
half the total bandwidth). The subsequent RG results, 
in particular the size of the BCDW region, depend on 
the choice of A cu toff and hence do not definitively an- 
swer the question whether the BCDW phase is intrinsic 
in EHM at half-filling. 

The virtue of the one-loop functional RG we develop 
and employ below is that it captures the BCDW phase 
in a systematic manner without ad hoc manipulations. 
The key point is that, while we truncate the flow equa- 
tions to order g 2 as in standard one-loop calculations, 
we maintain full momentum dependence of the interac- 
tion vertices. So instead of solving the RG flow equa- 
tions for four couplings g\, g 2 , 53 and 174, we write the 
functional RG equation for g(k\, k 2l ^3), where k±, A^, k 3 
can be anywhere in the Brillouin zone. Although the ra- 
dial dependence is formally irrelevant and the correpond- 
ing terms will eventually flow to zero, their effect may 
be finite when the energy cutoff is near the band edge, 
thereby breaking the accidental degeneracy. Other irrel- 
evant terms, such as higher-order vertices, are absent in 
the beginning of the flow and we neglect them altogether, 



just as in standard g-ology. We stress that our proce- 
dure should in general not qualitatively change the phase 
diagram-irrelevant operators will remain irrelevant-but 
may be crucial when an accidental degeneracy occurs. In 
this case very different phases may appear. 
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Figure 1: Discretization of the momenta in the Brillouin zone. 
This figure shows the case N = 10. 

Our functional RG equations for the one-dimensional 
EHM at one-loop follow closely the approach of Zanchi 
and Schulz 0], which itself is an adaptation of the 
two-dimensional RG for fermions 0] to the case of 
an arbitrary Fermi surface. This and other formula- 
tions of the functional RG have recently been applied 
to several two-dimensional interacting electron systems 
00 IH IH 0. The crucial difference is that we 
consider a finite number N of divisions of the magnitude 
of the momenta, while the two-dimensional calculations 
discretize the Fermi surface into angular patches. At one- 
loop level our equations become: 

dg(ki,k 2 ,k 3 ) 
dA 

/d 
dpj^[GA{p)GA(k)]g(ki,k2, k)g(p, k,k 3 ) 

f d 

- / dp—[G A (p)GA{qi)]g{p, k 2 ,qi)g(ki,qi,k 3 ) 
f d 

- dp—[GA(p)G A (q2)][-2g(k-L,p,q2)g(q2,k 2 ,k 3 ) 

+ g(p, h,Q2)g(Q2, k2,k 3 )+g(ki,p, q 2 )g{k 2 , q2,k 3 )], (2) 

where k — k\ + k 2 — p, qi = p + k 3 — k\ , q 2 = p + k 3 — 
k 2 , p = (p,uj), fdp = f dpYsu 1 /^ 271 ^)' and Ga is the 
propagator with cutoff A. 

The high-energy modes are integrated from the full 
bandwidth Ao (both positive and negative) to A, towards 
the Fermi surface. The cutoff A is parameterized by the 
RG parameter £ as A = Ao exp(— £). The initial condition 
for g{ki,k 2 ,k 3 ) is given by the Fourier transform of the 
U and V interaction terms. A direct analytical solution 
of the functional RG equation does not seem possible, 
and we use numerical calculations for solving the cou- 
pled integral-differential equations. For this purpose, the 
Brillouin zone is divided into N segments. Fig. ^ shows 
the discretization scheme for N = 10. 

We can follow the flows of the couplings to determine 
the dominant instability, but a more definitive answer 
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is given by comparing the susceptibilities corresponding 
to different broken symmetry states. We consider the 
susceptibilities of SDW, CDW, BSDW, and BCDW in 
the long- wavelength limit. Their general form is 

>,(3) 



where pi is the momentum at energy £j, J D(l,2) = 

i|f 1 |>A^l' / (^)i| S2 | > A C ^2J'(6)E < r 1 , CT2 ^i^) and J(0 
is the Jacobian for the coordinate transformation from fc 
to £ fc . For 5 = STW and (5 = BSTW:s T = 1, sj. = -1. 
For 8 = CDW and 8 = BCDW: s T = l,sj. = 1. For 
<5 = STW and 5 = CDW:f(p) = 1. For (5 = BSTW 
and 5 = BCDW: f(p) — sin(p). In momentum space, 
the difference between site and bond ordering is just in 
the form factor, which is s-wave for site orderings and 
p-wave for bond orderings. 
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Figure 2: The flows of SDW (solid lines), CDW (dotted lines), 
BSDW(dashed lines), BCDW(dotted-dashed lines) suscepti- 
bilities as function of t for U = 1, and V = 0.46,0.54,0.62. 
The left column is for N = 2, corresponding to standard g- 
ology, and the right column is for N = 50. 

Under the RG procedure, the susceptibilities also flow, 
with the dependence on A appearing both in the integra- 
tion and in the flow of the expectation value (...). The 



dominant instability is determined by the most divergent 
susceptibility as I is increased. The RG equations for the 
susceptibilities are, 



dA 
dA 



dp—[G A (p)G A (p)}(z{(p)) 2 , 



(4) 
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where p = p + (tt,0). For 8 = SDW and 6 = BSDW: 
g s {p',p) = ~g(p + n,p',p). For 8 = CDW and 8 = 
BCDW: g s (p',p) = 2g(p>,p + n,p)-g(p + 7r,p' 1 p). The 
function Z s (p) is the effective vertex in the definition 
for the susceptibility \ 5 ■ Its initial condition is I for 
SDW and CDW, and sin(p) for BSDW and BCDW. The 
RG equations for susceptibilities are solved with initial 
condition Xa=a i 71 ) = 0- 

It is instructive to compare the difference in the suscep- 
tibility flows for the (conventional) case of two scattering 
points N = 2 and for cases of multiple scattering points 
N > 2. We focus on the case U = 1 (in units of t) which 
is in the weak-coupling regime (U < Aq = 2t). In the left 
panels of Fig. we show the flows of susceptibilities for 
N = 2 i. e. for standard g-ology. In the right panels are 
results for functional RG with N = 50. For each case, re- 
sults are shown for V = 0.46, 0.54 and 0.62. These three 
values were chosen to cover the SDW, BCDW and CDW 
phases around the U = 2V line. 

First, we note that for N = 2, the susceptibilities tend 
to diverge more quickly with i. This is because, in effect, 
all the renormalization corrections to the irrelevant cou- 
plings (involving radial excursions away from the Fermi 
points) have been assigned to the marginal ones (gx, 
g4 of standard g-ology). This significantly enhances the 
rate of increase of the couplings and of the susceptibili- 
ties. Second, for N = 2 the SDW dominates for V = 0.46 
where U > 2V, and CDW dominates when U < 2V, as 
can be seen for V = 0.54 and 0.62. For N = 2 all the 
density wave susceptibilities are degenerate at U = 2V. 
Therefore, for N = 2 there is no finite region of BCDW 
phase. 

Importantly, for N > 2, we find that the BCDW 
suscetibility is dominant in a finite range around U = 2V. 
This is shown in the Fig. [21 for N — 50, for the case 
U = 1,V = 0.54. For smaller (larger) V the system is 
in the SDW (CDW) phase, as predicted by the standard 
g-ology. The pattern of SDW-BCDW-CDW for increas- 
ing V at a fixed U can be obtained by having only 6 
scattering points along the band. 

For our results to be reliable, it is important that 
they converge to a fixed result as N increases. Fig. 
shows that the phase boundaries converge quickly with 
N. Therefore, the N — 50 results presented in Fig. [21 
have reached the large N limit. 

The phase diagram we obtain is shown in Fig. 01 By 
using the momentum-dependent functional RG, we have 
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Figure 3: The shifts of the SDW-BCDW and BCDW-CDW 
phase boundaries with the number of momentum divisions N. 
One can note that the phase boundaries converge quickly with 
N, and they do not have any significant change for N > 30. 



well-controlled RG method at one loop level, with no ad- 
ditional interaction terms or ad hoc approximations. The 
width of this phase increases with U, and continues to ex- 
pand until the RG procedure breaks down. We have ver- 
ified that as N is increased the phase boundaries become 
independent of the number of divisions. Our results con- 
firm that a BCDW phase emerges spontaneously in the 
EHM at half-filling and clarify why this result has eluded 
earlier standard g-ology and bosonization techniques. 

We thank A. H. Castro Neto and A. Sandvik for useful 
discussions and Boston University for financial support. 



confirmed that the BCDW phase extends to very weak 
coupling (at least down to U ~ Ao/10) and expands with 
increasing U. This regime is difficult to access via QMC 
studies, and there has also been a controversy regarding 
the DMRG results in this limit. Intriguing phenomena in 
the strong-coupling regime, such as the shrinkage of the 
BCDW phase can not be reliably studied by this weak- 
coupling FRG method. This method takes into account 
irrelevant terms that couple spin and charge degrees of 
freedom, which has been argued in Ref.Q to be the cause 
for the disappearance of the BCDW phase at strong cou- 
pling. Nevertheless, the RG expansion is only valid for 
weak couplings and we focus on this limit in the present 
work. 
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Figure 4: The phase diagram of one- dimensional EHM at 
half- filling in the weak coupling regime. 

In summary, we have studied the phase diagram of the 
EHM at half-filling using a generalization of the stan- 
dard g-ology to a momentum-dependent, functional RG. 
In conventional terms, our approach includes formally 
irrelevant terms corresponding to interaction vertices in- 
volving electrons at high momenta. In the present case, 
our procedure changes the phase diagram qualitatively 
because it breaks an accidental symmetry between the 
backscattering (gi) and the Umklapp (53) processes at 
U = 2V. The full momentum dependence is included 
in a systematic way by discretizing the momenta in the 
Brillouin zone into N divisions. We obtain the BCDW 
phase near U = 2V by employing this consistent and 
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